precision highp float;varying vec2 vVertTexCoord;
uniform sampler2D source;
uniform sampler2D mask;
uniform vec4 addColor;
uniform vec3 random;
uniform bool active;

    #ifndef SPECTRAL
    #define SPECTRAL

    const int SPECTRAL_SIZE = 38;
    const float SPECTRAL_GAMMA = 2.4;
    const float SPECTRAL_EPSILON = 0.0001;

    float spectral_uncompand(float x) {
    return (x < 0.04045) ? x / 12.92 : pow((x + 0.055) / 1.055, SPECTRAL_GAMMA);
    }

    float spectral_compand(float x) {
    return (x < 0.0031308) ? x * 12.92 : 1.055 * pow(x, 1.0 / SPECTRAL_GAMMA) - 0.055;
    }

    vec3 spectral_srgb_to_linear(vec3 srgb) {
        return vec3(spectral_uncompand(srgb[0]), spectral_uncompand(srgb[1]), spectral_uncompand(srgb[2]));
    }

    vec3 spectral_linear_to_srgb(vec3 lrgb) {
        return clamp(vec3(spectral_compand(lrgb[0]), spectral_compand(lrgb[1]), spectral_compand(lrgb[2])), 0.0, 1.0);
    }

    void spectral_upsampling(vec3 lrgb, out float w, out float c, out float m, out float y, out float r, out float g, out float b) {
        w = min(lrgb.r, min(lrgb.g, lrgb.b));

        lrgb -= w;

        c = min(lrgb.g, lrgb.b);
        m = min(lrgb.r, lrgb.b);
        y = min(lrgb.r, lrgb.g);
        r = min(max(0., lrgb.r - lrgb.b), max(0., lrgb.r - lrgb.g));
        g = min(max(0., lrgb.g - lrgb.b), max(0., lrgb.g - lrgb.r));
        b = min(max(0., lrgb.b - lrgb.g), max(0., lrgb.b - lrgb.r));
    }

    void spectral_linear_to_reflectance(vec3 lrgb, inout float R[SPECTRAL_SIZE]) {
        float w, c, m, y, r, g, b;
        
        spectral_upsampling(lrgb, w, c, m, y, r, g, b);
        
        R[0] = max(SPECTRAL_EPSILON, w + c * 0.96853629 + m * 0.51567122 + y * 0.02055257 + r * 0.03147571 + g * 0.49108579 + b * 0.97901834);
        R[1] = max(SPECTRAL_EPSILON, w + c * 0.96855103 + m * 0.54015520 + y * 0.02059936 + r * 0.03146636 + g * 0.46944057 + b * 0.97901649);
        R[2] = max(SPECTRAL_EPSILON, w + c * 0.96859338 + m * 0.62645502 + y * 0.02062723 + r * 0.03140624 + g * 0.40165780 + b * 0.97901118);
        R[3] = max(SPECTRAL_EPSILON, w + c * 0.96877345 + m * 0.75595012 + y * 0.02073387 + r * 0.03119611 + g * 0.24490420 + b * 0.97892146);
        R[4] = max(SPECTRAL_EPSILON, w + c * 0.96942204 + m * 0.92826996 + y * 0.02114202 + r * 0.03053888 + g * 0.06826880 + b * 0.97858555);
        R[5] = max(SPECTRAL_EPSILON, w + c * 0.97143709 + m * 0.97223624 + y * 0.02233154 + r * 0.02856855 + g * 0.02732883 + b * 0.97743705);
        R[6] = max(SPECTRAL_EPSILON, w + c * 0.97541862 + m * 0.98616174 + y * 0.02556857 + r * 0.02459485 + g * 0.01360600 + b * 0.97428075);
        R[7] = max(SPECTRAL_EPSILON, w + c * 0.98074186 + m * 0.98955255 + y * 0.03330189 + r * 0.01929520 + g * 0.01000187 + b * 0.96663223);
        R[8] = max(SPECTRAL_EPSILON, w + c * 0.98580992 + m * 0.98676237 + y * 0.05185294 + r * 0.01423112 + g * 0.01284127 + b * 0.94822893);
        R[9] = max(SPECTRAL_EPSILON, w + c * 0.98971194 + m * 0.97312575 + y * 0.10087639 + r * 0.01033111 + g * 0.02636635 + b * 0.89937713);
        R[10] = max(SPECTRAL_EPSILON, w + c * 0.99238027 + m * 0.91944277 + y * 0.24000413 + r * 0.00765876 + g * 0.07058713 + b * 0.76070164);
        R[11] = max(SPECTRAL_EPSILON, w + c * 0.99409844 + m * 0.32564851 + y * 0.53589066 + r * 0.00593693 + g * 0.70421692 + b * 0.46420440);
        R[12] = max(SPECTRAL_EPSILON, w + c * 0.99517200 + m * 0.13820628 + y * 0.79874659 + r * 0.00485616 + g * 0.85473994 + b * 0.20123039);
        R[13] = max(SPECTRAL_EPSILON, w + c * 0.99576545 + m * 0.05015143 + y * 0.91186529 + r * 0.00426186 + g * 0.95081565 + b * 0.08808402);
        R[14] = max(SPECTRAL_EPSILON, w + c * 0.99593552 + m * 0.02912336 + y * 0.95399623 + r * 0.00409039 + g * 0.97170370 + b * 0.04592894);
        R[15] = max(SPECTRAL_EPSILON, w + c * 0.99564041 + m * 0.02421691 + y * 0.97137099 + r * 0.00438375 + g * 0.97651888 + b * 0.02860373);
        R[16] = max(SPECTRAL_EPSILON, w + c * 0.99464769 + m * 0.02660696 + y * 0.97939505 + r * 0.00537525 + g * 0.97429245 + b * 0.02060067);
        R[17] = max(SPECTRAL_EPSILON, w + c * 0.99229579 + m * 0.03407586 + y * 0.98345207 + r * 0.00772962 + g * 0.97012917 + b * 0.01656701);
        R[18] = max(SPECTRAL_EPSILON, w + c * 0.98638762 + m * 0.04835936 + y * 0.98553736 + r * 0.01366120 + g * 0.94258630 + b * 0.01451549);
        R[19] = max(SPECTRAL_EPSILON, w + c * 0.96829712 + m * 0.00011720 + y * 0.98648905 + r * 0.03181352 + g * 0.99989207 + b * 0.01357964);
        R[20] = max(SPECTRAL_EPSILON, w + c * 0.89228016 + m * 0.00008554 + y * 0.98674535 + r * 0.10791525 + g * 0.99989891 + b * 0.01331243);
        R[21] = max(SPECTRAL_EPSILON, w + c * 0.53740239 + m * 0.85267882 + y * 0.98657555 + r * 0.46249516 + g * 0.13823139 + b * 0.01347661);
        R[22] = max(SPECTRAL_EPSILON, w + c * 0.15360445 + m * 0.93188793 + y * 0.98611877 + r * 0.84604333 + g * 0.06968113 + b * 0.01387181);
        R[23] = max(SPECTRAL_EPSILON, w + c * 0.05705719 + m * 0.94810268 + y * 0.98559942 + r * 0.94275572 + g * 0.05628787 + b * 0.01435472);
        R[24] = max(SPECTRAL_EPSILON, w + c * 0.03126539 + m * 0.94200977 + y * 0.98507063 + r * 0.96860996 + g * 0.06111561 + b * 0.01479836);
        R[25] = max(SPECTRAL_EPSILON, w + c * 0.02205445 + m * 0.91478045 + y * 0.98460039 + r * 0.97783966 + g * 0.08987709 + b * 0.01515250);
        R[26] = max(SPECTRAL_EPSILON, w + c * 0.01802271 + m * 0.87065445 + y * 0.98425301 + r * 0.98187757 + g * 0.13656016 + b * 0.01540513);
        R[27] = max(SPECTRAL_EPSILON, w + c * 0.01613460 + m * 0.78827548 + y * 0.98403909 + r * 0.98377315 + g * 0.22169624 + b * 0.01557233);
        R[28] = max(SPECTRAL_EPSILON, w + c * 0.01520947 + m * 0.65738359 + y * 0.98388535 + r * 0.98470202 + g * 0.32176956 + b * 0.01565710);
        R[29] = max(SPECTRAL_EPSILON, w + c * 0.01475977 + m * 0.59909403 + y * 0.98376116 + r * 0.98515481 + g * 0.36157329 + b * 0.01571025);
        R[30] = max(SPECTRAL_EPSILON, w + c * 0.01454263 + m * 0.56817268 + y * 0.98368246 + r * 0.98537114 + g * 0.48361920 + b * 0.01571916);
        R[31] = max(SPECTRAL_EPSILON, w + c * 0.01444459 + m * 0.54031997 + y * 0.98365023 + r * 0.98546685 + g * 0.46488579 + b * 0.01572133);
        R[32] = max(SPECTRAL_EPSILON, w + c * 0.01439897 + m * 0.52110241 + y * 0.98361309 + r * 0.98550011 + g * 0.47440306 + b * 0.01572502);
        R[33] = max(SPECTRAL_EPSILON, w + c * 0.01437620 + m * 0.51041094 + y * 0.98357259 + r * 0.98551031 + g * 0.48576990 + b * 0.01571717);
        R[34] = max(SPECTRAL_EPSILON, w + c * 0.01436343 + m * 0.50526577 + y * 0.98353856 + r * 0.98550741 + g * 0.49267971 + b * 0.01571905);
        R[35] = max(SPECTRAL_EPSILON, w + c * 0.01435687 + m * 0.50255080 + y * 0.98351247 + r * 0.98551323 + g * 0.49625685 + b * 0.01571059);
        R[36] = max(SPECTRAL_EPSILON, w + c * 0.01435370 + m * 0.50126452 + y * 0.98350101 + r * 0.98551563 + g * 0.49807754 + b * 0.01569728);
        R[37] = max(SPECTRAL_EPSILON, w + c * 0.01435408 + m * 0.50083021 + y * 0.98350852 + r * 0.98551547 + g * 0.49889859 + b * 0.01570020);
    }

    vec3 spectral_xyz_to_srgb(vec3 xyz) {
        mat3 XYZ_RGB;

        XYZ_RGB[0] = vec3( 3.24306333, -1.53837619, -0.49893282);
        XYZ_RGB[1] = vec3(-0.96896309,  1.87542451,  0.04154303);
        XYZ_RGB[2] = vec3( 0.05568392, -0.20417438,  1.05799454);
        
        float r = dot(XYZ_RGB[0], xyz);
        float g = dot(XYZ_RGB[1], xyz);
        float b = dot(XYZ_RGB[2], xyz);

        return spectral_linear_to_srgb(vec3(r, g, b));
    }

    vec3 spectral_reflectance_to_xyz(float R[SPECTRAL_SIZE]) {
        vec3 xyz = vec3(0.0);
        
        xyz +=  R[0] * vec3(0.00006469, 0.00000184, 0.00030502);
        xyz +=  R[1] * vec3(0.00021941, 0.00000621, 0.00103681);
        xyz +=  R[2] * vec3(0.00112057, 0.00003101, 0.00531314);
        xyz +=  R[3] * vec3(0.00376661, 0.00010475, 0.01795439);
        xyz +=  R[4] * vec3(0.01188055, 0.00035364, 0.05707758);
        xyz +=  R[5] * vec3(0.02328644, 0.00095147, 0.11365162);
        xyz +=  R[6] * vec3(0.03455942, 0.00228226, 0.17335873);
        xyz +=  R[7] * vec3(0.03722379, 0.00420733, 0.19620658);
        xyz +=  R[8] * vec3(0.03241838, 0.00668880, 0.18608237);
        xyz +=  R[9] * vec3(0.02123321, 0.00988840, 0.13995048);
        xyz += R[10] * vec3(0.01049099, 0.01524945, 0.08917453);
        xyz += R[11] * vec3(0.00329584, 0.02141831, 0.04789621);
        xyz += R[12] * vec3(0.00050704, 0.03342293, 0.02814563);
        xyz += R[13] * vec3(0.00094867, 0.05131001, 0.01613766);
        xyz += R[14] * vec3(0.00627372, 0.07040208, 0.00775910);
        xyz += R[15] * vec3(0.01686462, 0.08783871, 0.00429615);
        xyz += R[16] * vec3(0.02868965, 0.09424905, 0.00200551);
        xyz += R[17] * vec3(0.04267481, 0.09795667, 0.00086147);
        xyz += R[18] * vec3(0.05625475, 0.09415219, 0.00036904);
        xyz += R[19] * vec3(0.06947040, 0.08678102, 0.00019143);
        xyz += R[20] * vec3(0.08305315, 0.07885653, 0.00014956);
        xyz += R[21] * vec3(0.08612610, 0.06352670, 0.00009231);
        xyz += R[22] * vec3(0.09046614, 0.05374142, 0.00006813);
        xyz += R[23] * vec3(0.08500387, 0.04264606, 0.00002883);
        xyz += R[24] * vec3(0.07090667, 0.03161735, 0.00001577);
        xyz += R[25] * vec3(0.05062889, 0.02088521, 0.00000394);
        xyz += R[26] * vec3(0.03547396, 0.01386011, 0.00000158);
        xyz += R[27] * vec3(0.02146821, 0.00810264, 0.00000000);
        xyz += R[28] * vec3(0.01251646, 0.00463010, 0.00000000);
        xyz += R[29] * vec3(0.00680458, 0.00249138, 0.00000000);
        xyz += R[30] * vec3(0.00346457, 0.00125930, 0.00000000);
        xyz += R[31] * vec3(0.00149761, 0.00054165, 0.00000000);
        xyz += R[32] * vec3(0.00076970, 0.00027795, 0.00000000);
        xyz += R[33] * vec3(0.00040737, 0.00014711, 0.00000000);
        xyz += R[34] * vec3(0.00016901, 0.00006103, 0.00000000);
        xyz += R[35] * vec3(0.00009522, 0.00003439, 0.00000000);
        xyz += R[36] * vec3(0.00004903, 0.00001771, 0.00000000);
        xyz += R[37] * vec3(0.00002000, 0.00000722, 0.00000000);

        return xyz;
    }

    float spectral_linear_to_concentration(float l1, float l2, float t) {
        float t1 = l1 * pow(1.0 - t, 2.0);
        float t2 = l2 * pow(t, 2.0);

        return t2 / (t1 + t2);
    }

    vec3 spectral_mix(vec3 color1, vec3 color2, float t) {
        vec3 lrgb1 = spectral_srgb_to_linear(color1);
        vec3 lrgb2 = spectral_srgb_to_linear(color2);

        float R1[SPECTRAL_SIZE];
        float R2[SPECTRAL_SIZE];

        spectral_linear_to_reflectance(lrgb1, R1);
        spectral_linear_to_reflectance(lrgb2, R2);

        float l1 = spectral_reflectance_to_xyz(R1)[1];
        float l2 = spectral_reflectance_to_xyz(R2)[1];

        t = spectral_linear_to_concentration(l1, l2, t);

        float R[SPECTRAL_SIZE];

        for (int i = 0; i < SPECTRAL_SIZE; i++) {
        float KS = (1.0 - t) * (pow(1.0 - R1[i], 2.0) / (2.0 * R1[i])) + t * (pow(1.0 - R2[i], 2.0) / (2.0 * R2[i]));
        float KM = 1.0 + KS - sqrt(pow(KS, 2.0) + 2.0 * KS);

        //Saunderson correction
        // let S = ((1.0 - K1) * (1.0 - K2) * KM) / (1.0 - K2 * KM);

        R[i] = KM;
        }

        return spectral_xyz_to_srgb(spectral_reflectance_to_xyz(R));
    }

    vec4 spectral_mix(vec4 color1, vec4 color2, float t) {
        return vec4(spectral_mix(color1.rgb, color2.rgb, t), mix(color1.a, color2.a, t));
    }

    #endif

// psrdnoise (c) Stefan Gustavson and Ian McEwan,
float psrdnoise(vec2 x, vec2 period, float alpha, out vec2 gradient) {
    vec2 uv = vec2(x.x+x.y*0.5, x.y);
    vec2 i0 = floor(uv), f0 = fract(uv);
    float cmp = step(f0.y, f0.x);
    vec2 o1 = vec2(cmp, 1.0-cmp);
    vec2 i1 = i0 + o1, i2 = i0 + 1.0;
    vec2 v0 = vec2(i0.x - i0.y*0.5, i0.y);
    vec2 v1 = vec2(v0.x + o1.x - o1.y*0.5, v0.y + o1.y);
    vec2 v2 = vec2(v0.x + 0.5, v0.y + 1.0);
    vec2 x0 = x - v0, x1 = x - v1, x2 = x - v2;
    vec3 iu, iv, xw, yw;
    if(any(greaterThan(period, vec2(0.0)))) {
        xw = vec3(v0.x, v1.x, v2.x);
        yw = vec3(v0.y, v1.y, v2.y);
        if(period.x > 0.0)
        xw = mod(vec3(v0.x, v1.x, v2.x), period.x);
        if(period.y > 0.0)
        yw = mod(vec3(v0.y, v1.y, v2.y), period.y);
        iu = floor(xw + 0.5*yw + 0.5); iv = floor(yw + 0.5);
    } else {
        iu = vec3(i0.x, i1.x, i2.x); iv = vec3(i0.y, i1.y, i2.y);
    }
    vec3 hash = mod(iu, 289.0);
    hash = mod((hash*51.0 + 2.0)*hash + iv, 289.0);
    hash = mod((hash*34.0 + 10.0)*hash, 289.0);
    vec3 psi = hash*0.07482 + alpha;
    vec3 gx = cos(psi); vec3 gy = sin(psi);
    vec2 g0 = vec2(gx.x, gy.x);
    vec2 g1 = vec2(gx.y, gy.y);
    vec2 g2 = vec2(gx.z, gy.z);
    vec3 w = 0.8 - vec3(dot(x0, x0), dot(x1, x1), dot(x2, x2));
    w = max(w, 0.0); vec3 w2 = w*w; vec3 w4 = w2*w2;
    vec3 gdotx = vec3(dot(g0, x0), dot(g1, x1), dot(g2, x2));
    float n = dot(w4, gdotx);
    vec3 w3 = w2*w; vec3 dw = -8.0*w3*gdotx;
    vec2 dn0 = w4.x*g0 + dw.x*x0;
    vec2 dn1 = w4.y*g1 + dw.y*x1;
    vec2 dn2 = w4.z*g2 + dw.z*x2;
    gradient = 10.9*(dn0 + dn1 + dn2);
    return 10.9*n;
}

vec4 generic_desaturate(vec3 color, float factor) {
    vec3 lum = vec3(0.299, 0.587, 0.114);
    vec3 gray = vec3(dot(lum, color));
    return vec4(mix(color, gray, factor), 1.0);
}

float rand(vec2 co, float a, float b, float c) {
    return fract(sin(dot(co, vec2(a, b))) * c);
}

void main() {
    vec4 maskColor = texture2D(mask, vVertTexCoord);
    if (maskColor.r > 0.0) {
        vec2 r=vec2(12.9898,78.233), a=vec2(7.9898,58.233), l=vec2(17.9898,3.233);
        float d=rand(vVertTexCoord,r.x,r.y,43358.5453)*2.-1., x=rand(vVertTexCoord,a.x,a.y,43213.5453)*2.-1., s=rand(vVertTexCoord,l.x,l.y,33358.5453)*2.-1.;

        const vec2 p = vec2(0.0);
        vec2 g;

        vec4 pigment;
        
        if (active) {
            float n = psrdnoise(vVertTexCoord * 5., p, 10.0 * random.x, g);
            float n2 = psrdnoise(vVertTexCoord * 5., p, 10.0 * random.y, g);
            float n3 = psrdnoise(vVertTexCoord * 5., p, 10.0 * random.z, g);
            float tt = abs(addColor.x - addColor.y - addColor.z);
            float n4 = 0.25 + 0.25 * psrdnoise(vVertTexCoord * 4., p, 3.0 * random.x, g);
            pigment = vec4(generic_desaturate(addColor.xyz,n4).xyz + vec3(n,n2,n3) * 0.03 * tt, 1.0);
        } 
        else {
            pigment = vec4(addColor.xyz,1.0);
        }

        float darken_above = 0.7;
        if (maskColor.a > darken_above) {
            float blacken = 0.5 * (maskColor.a - darken_above);
            pigment = pigment * (1. - blacken) - vec4(0.5) * blacken;
        }

        vec3 mixedColor = spectral_mix(texture2D(source,vVertTexCoord).xyz, pigment.xyz, 0.9 * maskColor.a);
        gl_FragColor = vec4(mixedColor + 0.01*vec3(d,x,s),1.);
    }
}